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We consider a single particle which is bound by a central potential and obeys the Dirac equation 
in d dimensions. We first apply the asymptotic iteration method to recover the known exact solutions 
for the pure Coulomb case. For a screened-Coulomb potential and for a Coulomb plus linear potential 



with linear scalar confinement, the method is used to obtain accurate approximate solutions for both 
eigenvalues and wave functions. 
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1. Introduction 

The Dirac equation plays a fundamental role in relativistic quantum mechanics. The equation 
can be solved exactly for a very few potentials. An early and very detailed analysis of the Dirac 
spectrum for central potentials has been given by Rose and Newton [1]. The solution for the pure 
Coulomb field is well-known. Exact solutions are known for some other specific cases, such as the 
Woods-Saxon potential [2]. It is also possible to find classes of exact solutions under very special con- 
ditions [3-7]. Unlike the corresponding Schrodinger operator, the Dirac Hamiltonian is not bounded 
below, and its spectrum can not be defined variationally. In spite of this it is still possible to find 
ways to use the variational idea, for example by a saddle-point analysis [8] and by an analysis of 
the Schrodinger limit in a oscillator basis [9]. Some progress has been made in the establishment of 
comparison theorems for the Dirac equation without invoking a variational assumption at all [10,11], 
and this has allowed spectral envelope methods to be used; however, at present, such results apply 
only to node-free states. Theorems of the Levinson type have now been proved for the Dirac equa- 
tion in two, three, and d dimensions [12-18]. Solutions of the Dirac equation with shape-invariant 
potentials have been found by a variety of exact and approximate methods [19-29]. Meanwhile, nu- 
merical methods for solving the Dirac equation are continually sought, such as a recent approach 
by a mapped Fourier grid method [30]. Effective non- variational approximation methods therefore 
remain an important area for investigation. 

In this paper we first consider the Dirac equation for central potentials in d dimensions. This 
problem was formulated some decades ago, by Joseph [31]. More recently, Yu Jiang [32], for example, 
has studied the problem and has obtained a pair of radial equations similar to the well-known 
case of three dimensions; these equations can be solved exactly for the pure Coulomb field. The 
d -dimensional angular momentum problem for the Dirac equation has also been studied by group- 
theoretical methods [33]. We use these results in section 2. The main purpose of the present paper 
is to apply the asymptotic iteration method (AIM) to these central-field problems. AIM was first 
developed by the present authors [34] for solving second order linear differential equations, including 
Schrodinger's. Later the method has been applied to a variety of problems [35-37]. In preparation 
for the Dirac application, we first extend the method in section 3 to treat systems of homogeneous 
linear differential equations. In section 4 the case of the Dirac Coulomb problem in d dimensions is 
then treated and solved exactly. In sections 5 and 6 we study screened Coulomb problems, and also 
the linear plus Coulomb potential. For the latter problem, the scalar linear part must dominate the 
vector linear part in order for discrete eigenvalues to exist [38, 39]. 
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2. The Dirac equation in d dimensions 

The Dirac equation for a central field in d dimensions can be written with h = c = 1 as 

= HV, H = J2 <XjPj + P(m + U{r)) + V(r), (2.1) 

where m is the mass of the particle, V(r) is a spherically symmetric vector potential, U(r) is 
a spherically symmetric scalar potential, and {ctj} and (3 are the usual Dirac matrices satisfying 
anti-commutation relations (the identity matrix is implied after the vector potential term V(r) ). 
After some algebraic calculations (details can be found for example in [32]), one obtains the following 
first-order linear coupled differential equations 

^ = -^G + (E + m- V(r) + U(r))F (2.2) 

4— = —{E — m — V(r) - U(r))G + —F (2.3) 
dr r 

These are known as the radial Dirac equations in d dimensions, where kd = r(j + ^p) , and 
r = ±1. We note that the variable t is sometimes written ui , as, for example, in the book by 
Messiah [40] . As an example, we suppose that the particle is moving in a pure vector Coulomb field, 
that is to say, V(r) = — 7- and U(r) = . In this case Eqs.(2.2) and (2.3) can be written as follows 

dG =- k -±G+(m + E+*)F (2.4) 



dr r \ r 

d JL=( m -E-*)G + k -±F. (2.5) 
dr \ r ) r 

First of all, we have to obtain some asymptotic forms for G{r) and F(r) functions. At small r , 
G{r) can be written as the following Euler equation 

« !« + (fi^£) c (2 . 6) 

dr r dr \ r z J 

It is clear that the solution of Eq.(2.6) is: G(r) — r 7 , where 7 = \Jk^ — A 2 . When a similar 
analysis at small r is made for F(r) , one finds the same general results. At large r , we obtain 
the following asymptotic differential equation for both G(r) and F(r) : d = (in 2 — E 2 )H(r) , 
where H(r) is either G(r) or F(r). We conclude that H(r) ~ exp(— \Jm 2 — E 2 r) . We therefore 
adopt the following representations for these radial functions 

G(r) =Vm + E r 1 exp {-r\J m 2 - E 2 ^ (0i + <j> 2 ) (2.7) 
F{r) = ^m-E r 7 exp {-r^/m 2 - E 2 } (fa - <f> 2 ) (2.8) 
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We now substitute these equations into the Eqs.(2.4) and (2.5) and use the notation r = np to 
find: 

#1 ^l-i±lU-f^U (2.9) 



dp V P 



am 

dp 



b-k d \ fa -7 s 



(2.10) 



where r\ = 2 ^ J J 2 _ E2 > a = 2EAri and b = 2mAr\ . It is possible, of course, to solve these equations 
by using a power series method. However, we shall solve the equations by means of AIM, which was 
developed originally for second order linear differential equations. In the next section, we extend the 
scope of AIM to apply to first order linear coupled differential equations generally; then in section 4 
we apply the results obtained to the specific Dirac radial equations (2.9) and (2.10) above. 

3. Asymptotic iteration method for first order linear coupled differential equations 

We consider the following first order linear coupled differential equations 

01 = \o(x)(f>i +s (x)<j) 2 (3.1) 

02 = uo(x)<(>i +Po(x)<j)2, (3.2) 

where {/} represents the derivative with respect to x, and Ao(a;) , sq(x) , u)q(x), and po(x) arc 
sufficiently differentiablc in appropriate domains. If we differentiate (3.1) and (3.2) with respect to 
x , we find that 

4>'{ = \ 1 {x)<t> l +s 1 {x)4>2 (3.3) 
^'=wi(a:)0i+pi(a;)^2, (3.4) 

where 

Ai = Aq + Aq + sqloq, 
si = s' + XqSq + s po, 

Pi = Po+Pq + s o^o- 
Similarly, if we calculate the (n + 2) th derivative, n = 1,2,..., we have 

4>[ n+2) = A n+ i(a;)0i + s n+1 (x)(f> 2 (3.5) 
2 ™ +2) = w n+1 {x)(pi +p n+ i{x)(f>2, (3.6) 
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where 

Ki+i = A n + A„Aq + s n LUQ, 
s n +i = s' n + X n s + s n p , 

LU n+1 = U)' n + CJ„A + Pn^O, 
Pn+1 =Pn+V n S0+PnP0- 

From the ratio of the (n + 2) th and (n + l) th derivatives of <j>\ we have 



dx 



In 



(3.7) 



An exactly similar result can be obtained for <fi 2 ■ However, to solve the system given in (3.1) and 
(3.2), one of these conditions is sufficient. We now introduce the 'asymptotic' aspect of the method. 
If we have, for sufficiently large n , 



= — :=a, n=l,2,3,. 

A n +1 A n 



then (3.7) reduces to £ ln((^" +1) ) = ) w hich yields 



(3i 



X I X \ 

0<" +1) (x) = & exp ( | ^|^^) = CiA„ exp J {au a + X )dt , (3.9) 

where C\ is the integration constant. After substituting Eq.(3.9) into </>i™ +1 ' = \ n (x)<f>i + s n (x)(j)2 , 
we get 

4>i + a(x)4> 2 = C\ exp ^ y (aw + \ )dt 
Using Eqs.(3.10) and (3.2), we can obtain the general solution of <p2(x) as follows 

<fo{x) = exp iJ{po~ uJ ot)dt I C 2 + Ci J I ljo exp I J (A - p + 2w aO dT J cfe 



(3.10) 



(3.11) 



Once we have obtained (j> 2 {x) , it is easy to find <pi{x) by using one of the coupled equations, or 
directly from Eq.(3.10). 



4. Solution of the Dirac Coulomb problem 

We now turn back to our first principal application. If we compare Eqs.(2.9) and (2.10) with Eqs.(3.1) 
and (3.2), we see that X (p) = 1 - 2±r , s (p) = , wo(p) = h -=^ and p (p) = ^2 . By 

using our iteration formulas and the iteration termination condition given in Eq.(3.8), we find that 
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a = 7, 1 + 7, 2 + 7, ... this means that a = n + 7 , where n = 0, 1, 2, ... . In this case b satisfies the 
following relations 

n = 0; a = 7; 6 = — fc^ 



1; a=l+7; & = -fc d , ±Wl + 27 + fc2 



n 



= 2; a = 2 + 7; 6 = ±^3 + 27 + ^, ±^4 + 47 + ^ 

-3; a = 3 + 7; 6 = ±^5 + 2 7 + k 2 d , ±^8 + 4-/ + k 2 d , ±^9 + 6j + kj. 



In general we have, for a 



7, that b = ±\Jk\ + s(2n - s) + 2sj , where s = 0,1,2,. 



We 



know from Eqs.(2.9) and (2.10) that a — 2EAr\ , b = 2mAn, and r\ — 2 ^ r J 2 E ± 
these equations, we find the following two different expressions for the energy 



E = ± 



m 



. When we use 



(4.1) 



and 



E = ±mjl 



A 2 



k\ + s(2n - s) + 2s-/' 



(4.2) 



But these expressions must be equal. From their equality, s must be n or n + 2j . We know that 
both n and s are integers, but 7 is not an integer, so we find that s = n . In this case 6 becomes 
b = ±^Jk\ + n 2 + 2n7 , but when n = , b must be — kd ■ Thus the energy has the following form 

-1/2 



E = ±m 



1 



n+^k 2 - A 2 



(4.3) 



where we use 7 = y/k% - A 2 and n = 0, 1, 2, 3, ... . If we define the principal quantum number as 

2 



n r = n+ \kd\ — = 1,2,3,..., we recover the following well-known formula for the Coulomb energy 



E = ±m 



1 + 



n r - \k d \ + + ^k 2 - A 2 



-1/2 



(4.4) 



5. Coulomb wave functions 

In this section we obtain the Coulomb cigenfunctions by using Eqs.(3.14) and (3.11). Eq.(3.11) 
includes two independent solutions for <p2 {x) . The first factor of the expression in Eq.(3.11) generally 
represents the physical solution, so we use this factor as the wave-function generator for our model. 
If we rewrite this, we have 

<fo{x) = C 2 exp / (p - uj a)dt , (5.1) 
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where C2 is an integration constant which can be determined by normalization. If we use our 
iteration procedure, we find the following results for <p 2 ■ 

h(p) = 1, n = 

6 2 (p) = -(2 7 +l)(l-^- T ) ) n = l 

H p) = (27 + D(27 + 2)(1 - ^ + (27+1 ^ a7 + 2) ), " - ^ 

MP) = -(27 + D(27 + 2)(2 7 + 3)(1 - ^ + ^ - ^ - 2) - ^ - ^ ^ - g) ), n = 3. 

We see from these results that the general formula for <f>2 (p) can be written as follows 

2 (P) = (-1)" (2 ^ ) " )! C 2 iFi(-n,2 7 +l,p), (5.2) 

where 1-F1 is the confluent hypergeometric function; since the first argument is a negative integer, 
the function is a polynomial of degree n . We can now calculate 4>i (p) . For this task we use Eq. (3.10) , 
but like 4>2{p) , the solution of Eq.(3.10) has two parts, one of them is a polynomial and the other 
one is an infinite series: we choose the polynomial solution. Thus we find that 4>i(p) = —a(p)(j)2{p) , 
where a(p) = jf*- . If we calculate 4>i{p), using the above equation, we find the following results: 

MP) = 0, n = 
4>i(p) = (kd + b), n=l 

MP) = -(k d + 6)(27 + 1)(1 - 2^-), n = 2 

MP) - ft, + 6)(27 + D(27 + 2)(1 - ^ + (27+1 ^ 27 + 2) ), « - 3. 



We conclude from these results that 

Mp) = {-i) n+ \b + k d ) {2l+ { ^ l) - C2 ,r,(:i -,,.2- + i.//). (.",;; 



where 6 = ^j^ A E t ■ After obtaining (f>\ and <p2 , we can recover the radial functions G{r) and 
F(r) in the well-known form by using Eqs.(2.7) and (2.8). Thus, we have the complete solution of 
the Dirac equation for the Coulomb problem. In the next sections, we will discuss the solution of the 
Dirac equation for a screened-Coulomb potential and Coulomb plus linear potential with a linear 
scalar confinement. 

6. Dirac equation for a Screened-Coulomb potential in 3-dimensions 

In this section we will study a screened-Coulomb potential in 3-dimcnsions. We use the Mehta 
and Patil potential [41] which is suitable for large atoms. This potential is defined by 

V(r) = h — , 6.1 

r 1 + Xr 
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where Vi — Za , v 2 = (Z — l)a and A = O.^aZ 1 ' 3 , and a w 1/137.036 is the fine-structure 
constant. For this case, the radial Dirac equations read 

^- = --G+(m + E+^--W(r)]F (6.2) 
dr r \ r / 

^= ( m -E-^ + W(r))G+-F, (6.3) 
dr V r J r 

where W(r) = i+^ r 1 asymptotic behaviours of -F(r) and G(r) are the same as for the pure 
Coulombic potential. For this reason, G(r) and F(r) can be written as follows 

G(r) = r 7 cxp (-ry/ 'm 2 - E 2 } (fa + 2 ) (6.4) 
F(r) = r 7 cxp {-ryj m 2 - E" 1 } (fa - fa), (6.5) 

where 7 = yjk" 1 — v\ . After substituting these forms into into Eqs.(6.2) and (6.3), we have 

^ = ( m + ,-l) 4 , 1 -(E-W(r) + ^) <h (6.9) 



2 



E-W(r) + ^ ^\fa+(a-m-^fa, (6.7) 



v\ — k 

dr \ w p 

where cr = \/m 2 — E 2 . Now we can use the AIM procedure to obtain the eigenvalues. In preparation 
for the iteration process it helps first to remove the square-root expression in a . We do this by the 
following subsitutions 



m-E 2me „ / e 1 - 1 \ ,„ „. 

, cr = ^, E = m[- n . (6.8 

m + E 1 1 + e 2 ' V £ + 1 / 

If we use our iteration formulas and the iteration termination condition given in Eq.(3.8), we can 
construct the following equation, which corresponds to Eq.(3.8): 

S(r, e) = A„+is„ - s n+1 \ n = 0. (6.9) 

If the problem is exactly solvable, then 5(r, e) = 5(e) is independent of r and its vanishing gives 
us the exact results, as with the pure Coulomb problem discussed in section 4. In cases that are not 
exactly solvable in closed form, S(r, e) depends on both r and on e. We then solve the equation 
S(ro,e) = for a suitable fixed r = r point, which choice affects the convergence rate of the 
iteration (this choice is discussed in more detail in the next section). For the problem at hand we 
chose only one value, rg = 2, which fixed choice led to fast convergence in all cases. As the iteration 
number increases, the eigenvalue estimates become more accurate. The results for the ground state 
energies (i.e for ls 1 / 2 ) with m = 1 and various atomic numbers Z are presented in Table 1. These 
agree with the results obtained earlier in Rcf. [10]. 
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7. Dirac equation for the Coulomb plus linear potential with a linear scalar term 

As a further test of the method, we turn in this section to a problem quite different from that 
of atomic physics: we study the Dirac equation in 3-dimcnsions in the case that the vectorial part 
of the potential is Coulomb plus linear, and the scalar part is linear. We should like to point out 
that exactly similar calculations can be carried out in arbitrary dimensions d . For this problem the 
radial Dirac equation is written as follows 

[a.p + P(m + U(r)) + V(r)]^ = E^, (7.1) 

where V(r) = — ^ + B\r and U(r) = B 2 r , A is positive, and B 2 > B\ . The correspoding radial 
Dirac equations become 

^ = --G+(m + E+- + (B 2 -B 1 )r)F (7.2) 
dr r \ r J 

dF {A \ k 

— = (m-E-- + (B 2 + B 1 )r )G+-F, (7.3) 
dr \ r J r 

where k = r(j + j) and r = ±1. Before starting to calculate the energy eigenvalues, we have to 
consider the asymptotic behaviour of F(r) and G(r) at the boundaries. First we consider small 
r . Since the Coulomb potential dominates in this region, the asymptotic behaviour is the same as 
for the Hydrogenic problem: that is to say, F{r) and G{r) behave as r 7 , where 7 = \Jk 2 — A 2 . 
At large r , we find that H"(r) ~ (2(mB 2 + EB x )r + [B\ - B 2 )r 2 ) H{r) , where H{r) is G(r) or 
F(r) . Thus, H(r) - cxp (-ar - \(3r 2 ) , where (3 = B 2 - B 2 and a = mB ^ EB ^ . Thus we now 
write G(r) and F(r) in the following way 

G(r) = r 7 exp (-ar - ^(3r 2 ^j (fa + fa) (7.4) 

F(r) = r 7 exp (-ar - i/3r 2 ^ (fa - fa). (7.5) 

After substituting Eqs.(7.4) and (7.5) into (7.2) and (7.3), we find 



dr 
dfa 
dr 



-=(m + a+(l3 + B 2 )r - J) fa - (e + - B.r^j 



(7.6) 



2 ( E +±± 



B 1 r^jfa + ^a-m-^ + (f3-B 2 )r^fa. (7.7) 



When we compare these equations with (3.1) and (3.2) we see that A (r) = m + a + (rj + B 2 )r — ^ , 
s (r) = -E - + B x r , w (r) = E + - B x r and p (r) = a-m-^ + (r]- B 2 )r . We now 
have to solve S(E, r) = at a suitable ro point. Thus ro is a parameter of the method. For the 
present problem, we have found that any ro satisfying 1 < ro < 3 is satisfactory in the sense that 
the iteration sequence converges rapidly. For non-relativistic problems, we have found earlier that 
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r can be chosen as the peak of a simple scale-optimized trial function; a value of r found in this 
way for the the ground state is also effective for the excited states. However, since we don't have a 
suitable variational principle for the Dirac case, the 'scale' as represented by r , is chosen by the 
convergence criterion. For the specific example discussed here, we adopted the fixed value ro = 1.5 
and we obtained E for A = 1/2 , B 2 — 0.2, and B\ = 0.1 . The results are exhibited in Table 2 
(the spectroscopic labelling is explained below). When B 2 is much bigger than B\, AIM gives us 
more accurate results for small iteration numbers, but here, we calculate the energy eigenvalues for 
rather close values of B 2 and B\ , in order to test the effectiveness of the method. 

For the Dirac equation with central potentials I is not a good quantum number. However, 
a spectroscopic description of the states is still possible if we adopt the following convention [10]. 
We recall from section 2 that r = ±1. Meanwhile, the lower index I = 0, 1, 2, . . . of the spherical 
harmonic Y™ appearing [40] in the upper two components of the Dirac 4-spinor is related to j 
by t = j + |t, and the parity of the state is given by the formula P = (— iy + i T = (— l) 1 . If the 
number n = 1, 2, 3, . . . counts the eigenvalues for each given value of k = r(j + 1), we may then 
define the 'principle quantum number' for all central potentials as n r = n + I. The number I can 
then be represented by the usual atomic symbol {s,p, d, /..} . With this convention we label a state 
by n r Dj , where D = s,p, d, /..; In the non-relativistic large- m limit, this notation agrees with the 
usual Schrodinger description. 

8. Wavefunctions for the linear plus Coulomb problem 

It is possible to obtain approximate wavefunctions for the potential with Coulombic vectorial 
and scalar linear parts. We calculate the wavefunctions approximately by using AIM in the following 
way. As in section 4, we use the first part of Eq.(3.11) to generate the wavefunction. In this case, the 
function a(r)wo(r) which appears in the wave function generator can not be integrated analytically 
at every iteration, so instead of doing this, we first expand this function near r — and then integrate 
the representation to give an approximation for (j)2{r). Once we have <j}i{r) , it is straightforward 
to obtain (f>i(r) by using (j>i{r) — -a(r)^(r) . After obtaining <pi{r) and (fe{r) , we can find G(r) 
and F(r) by using Eqs.(7.4) and (7.5). Below we give an example of these calculations. We choose 
the 3p 3 / 2 state in Table 3 and find the following corresponding wave functions: 



\fe=o / 

where and bk coefficients are given in Table 2. We tabulate more coefficients than are needed, 
in order to demonstrate the stability of the method. In Figure 1 we show Cartesian plots of the 




(8.1) 



(8.2) 
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radial functions G(r) and F(r), and in Figure 2 we exhibit the Dirac spinor orbit [38] defined 
by (G(r) , F(r) ), r > 0. The availability of the wave-function approximations makes these tasks 
straightforward. 

9. Conclusion 

In this paper we have shown how AIM can be used to solve systems of two first order linear 
differential equations. In cases where the system represents an eigenvalue problem, the method yields 
the eigenf unctions and the eigenvalues. If the exact wave function may be factored in the form of an 
asymptotic wave function multiplied by a polynomial, the problem can be solved exactly. In other 
cases, an approximate solution is found by forcing the vanishing of a certain function 5(r,E) after 
a finite number of iterations at a fixed expansion point r = r . The range of values of ro which all 
lead to fast convergence is not narrow: for the problems discussed in this paper the range 1 < r < 3 
was satisfactory; for highly excited states with Coulomb-like potentials, we have found larger values 
to be better. 

In this paper we report applications of the method to bound states of the Dirac equation. First 
of all, as a test, the known exact solutions of the Coulomb problem in d dimensions were recovered. 
The method was then applied to find the spectrum and wave functions for a screcncd-Coulomb 
potential, and also for a very different problem, namely a linear plus Coulomb potential with a 
scalar linear confining term. In all cases the method yielded fast convergence to accurate solutions. 
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Table 1 Ground state eigenvalues E for the state k = — 1, j = 1/2 (with spectral description 
lsi ) for the screened Coulomb potential. The energies (E — l)m e , where m e = 511.004 KeV, are 
shown, along with corresponding accurate numerical values for comparison. 



z 






l\Ti lmorirfi 1 


20 


991560 


-4.3129 


-4.3157 


30 


0.979852 


-10.2957 


-10.2960 


40 


0.962675 


-19.0732 


-19.0732 


50 


0.939619 


-30.8549 


-30.8543 


60 


0.910139 


-45.9193 


-45.9189 


70 


0.873475 


-64.6548 


-64.6545 


80 


0.828543 


-87.6152 


-87.6148 
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Table 2 Eigenvalues E of the Dirac Hamiltonian for the linear plus Coulomb vector potential 
V(r) = —A/r + B\r and a linear scalar potential U(r) = B2 r in three dimensions, where A = 

0. 5, B2 = 0.2, B\ = 0.1 The results are given in dimensionlcss units corresponding to m = 

1. Accurate numerical results E num are shown for comparison: this accuracy was obtained with 20 
iterations. 



k 


n 


States 


E 




-1 







1.25819 


1.25819 


-1 


1 


2S1/2 


1.87575 


1.87575 


-1 


2 


3si/ 2 


2.29722 


2.29722 


1 





lPl/2 


1.70367 


1.70367 


1 


1 


2pi/ 2 


2.15272 


2.15272 


1 


2 


3^1/2 


2.51020 


2.51029 


-2 





2^3/2 


1.74683 


1.74683 


-2 


1 


3^3/2 


2.19096 


2.19096 


-2 


2 


4p 3 /2 


2.54480 


2.54486 


2 





24/2 


2.03889 


2.03889 


2 


1 


34/2 


2.41193 


2.41193 


2 


2 


44/2 


2.72766 


2.72762 


-3 





34/2 


2.04506 


2.04506 


-3 


1 


44/2 


2.42019 


2.42019 


-3 


2 


54/2 


2.73666 


2.73665 
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Table 3 The Coefficients a k and b k in Eqs(8.1) and (8.2) 



k 


ak 


bk 





1.7746 


-0.22540 


1 


3.34842 


-0.87777 


2 


2.58401 


-1.16405 


3 


0.89712 


-0.79148 


4 


-4.42384. 10" 2 


-0.30087 


5 


-0.20572 


-4.92821. 10" 2 


6 


-0.11569 


1. 18971. 10" 2 


7 


-3. 99651. 10~ 2 


1.06347.10" 2 


8 


-1.00424.10- 2 


3. 82882. 10" 3 


9 


-1.95538.10" 3 


9. 31883. 10- 4 


10 


-3. 04816. 10~ 4 


1. 71231. 10" 4 


11 


-3. 87895. 10~ 5 


2.48733. 10" 5 


12 


-4.07966. 10- 6 


2. 92628. 10" 6 


13 


-3. 57403. 10~ 7 


2. 82781. 10" 7 


14 


-2. 61991. 10~ 8 


2. 26316. 10" 8 


15 


-1. 61008. 10~ 9 


1. 50639. 1Q- 9 
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Figure 2. Dirac spinor orbit (F(r),G(r)), r > 0, in dimensionless units for the same example as 
in Figure 1 



